Fitting the dielectric response of collisionless plasmas by continued fractions 
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We present an approximation scheme for the dielectric response of thermal collisionless plasmas 
at arbitrary degeneracy. A T-fraction representation is obtained from the known expansions of 
the real part of the dielectric function for small and large arguments. The partial numerators and 
denominators of the continued fraction are generated by a modified Q-D algorithm. For several 
typical values of the degeneracy parameter 9, extensive tables for the expansion coefficients and the 
partial numerators and denominations are given allowing for an easy implementation of the fitting 
function. Also, an error analysis is performed. 



I. INTRODUCTION 

The dielectric response of a collisionless plasma is an ubiquitous quantity in plasma physics and in the many-body 
theory of Coulomb systems in general. For a fermionic system, it was first derived by Lindhard [l[ and is closely 
related to the dielectric function e(k, w) in random phase approximation (2j. In this way, the dielectric response is 
connected to collective effects such as screening, plasmons, and Landau damping. Its knowledge is also important for 
non-ideal plasmas, since the effects of non-ideality are typically parameterized by dynamic local field corrections with 
respect to the ideal response @, H| . 

A detailed study of the ideal dielectric response for fermionic systems has been reported by Arista and Brandt [![■ 
Here, we retain this restriction to fermions, although a treatment for bosonic systems is possible in a similar manner. 
In particular, we consider an one-component system with temperature T and density n interacting by the Coulomb 
potential. Introducing the function g(x) as 



ydy 



>o exp{y 2 /e~r]) 
the real part of the dielectric function for fermions is given by 
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with u = u)/(kvp), z = k/(2kp), Xq = (^fOb) \ and the degeneracy parameter 9 = ksT/Ef. Here, Ep,VF,kF 



1/3 



are the Fermi energy Ep = fi k F /(2m), Fermi velocity vf = hkp/m, and Fermi wave number kp = (3n n) 
respectively. The Bohr radius is denoed by as, f] the chemical potential obtained from 2/3 = -^1/2(7?) # 3 ^ 2 with the 
Fermi function The imaginary part of the dielectric function is known analytically for all degeneracies 
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For special cases, see the tables in Rcf. Q. 

The dielectric response enters a number of important physical observables such as optical properties 0] , the single- 
particle self-energy in GW approximation Q , and the dynamical collision frequency Q . A fast and reliable computa- 
tion is therefore of crucial importance to obtain these quantities. However, a direct evaluation of the integral Eq. ([1]) 
is often too slow in many applications. Thus, approximative analytical expressions for g{x) are desirable. 

Continued fractions appear in a variety of applications in theoretical physics 0. Also, they are of importance 
in approximation theory and are closely related to Pade approximants which in turn enjoy a number of interesting 
applications in various fields [1, In particular, there is a close correspondence between T- fractions and two 
point Pade approximations, see Refs. [Il| . Since an expansion of g(x) at z = and at z = 00 is known, the T-fractions 
representation suggests itself as a powerful approximation. Also, due to the compact form of the continued fraction, 
only a few fitting parameters have to be given. 
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II. APPROXIMATION SCHEME 

Following Arista and Brandt we can consider the expansion of g(x) for small and large values of x and obtain 

g(x) = 2H 1 (6)x + Ih 2 (9)x 3 + ^H^x 5 + ... -i-fr^x*" 1 + .. (4) 
3 5 Zi — 1 

9{X) = 3~x + d + 5x5 +■■■ + £j + ■•• (5) 

where is the Fermi integral of order v. i.e. 

/*OC 

FM = I dt YT - r , (6) 

and i being an integer, I an odd number. Note, that Hi are hyper-singular integrals which have to be regularized. 
Details can be found in appendix [Al 

Furthermore, the function g(x) is analytically known in the non-degenerate limit 9> 1 

9 {x) = 2 -e-^v(x/e^) , (7) 

and in the highly degenerate limit 6 <C 1 

1 + x 



g(x) = x + i (l - .t 2 ) In 



(8) 



V{x) is the so-called plasma dispersion function [ll 



1 f 00 e~ y 

V(x) = -=P / dy , (9) 

and closely related to the complex error function [13l|. P J indicates a Cauchy principal value integration. 

It has been shown by McCabe and Murphy [111, that there exists a close connection between a two-point Pade 
approximant and a T-fraction [Icj given in general by 

Fiz I F 2 z F 3 z , . 

with Fj, Gj being some constants and z being a complex variable. Here, we follow the notation introduced in Ref. [Til ]. 
In particular, once the expansions for z = and z = 00 are known to be of the form 



(11) 



z z 2 z 3 z fe+1 
and 

- [i- 2 z - H-3Z 2 - ... - H-kZ^ 1 + ... , (12) 

with expansion coefficient ^ and a complex variable z, the T-fraction representation can be generated by a Q- 
D algorithm provided that all expansion coefficients are different from zero. A modified Q-D algorithm has been 
developed by de Andrade et al. [la] to allow for zero expansion coefficients. 

Now, the above given expressions for g{x) are almost of the form necessary for applying the T-fraction scheme. We 
have to take care of the infinitesimal small imaginary contribution to the expansion at 00 

Imflr(a;) = |01n(l + cxp (r? - x 2 /6»)) . (13) 

by defining g{x) = g(x) — ilmg(x) and calculating the continued fraction representation for g(x). After generating 
the representation, we can go back to the original g(x). 
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Hi 


H 2 


H a 


H 4 


H 5 


0.1 


9.914(-1) 


-1.0186 


-0.3649 


-0.2482 


-0.2108 


0.2 


9.611(-1) 


-1.0892 


-0.4569 


-0.2997 


-0.1499 


0.3 


9.081(-1) 


-1.1533 


-0.4130 


-0.0950 


0.1268 


0.5 


7.790(-l) 


-1.1037 


-0.1000 


0.1734 


0.1217 


0.8 


6.116(-1) 


-0.8336 


0.1250 


0.0939 


-0.0075 


1.0 


5.289(-l) 


-0.6692 


0.1424 


0.0403 


-0.0133 


1.5 


3.888(-l) 


-0.3970 


0.0973 


-0.0005 


-0.0038 


2.0 


3.048(-l) 


-0.2552 


0.0581 


-0.0038 


-0.0007 



TABLE I: Expansion coefficients for the small x expansion of the function g(x), cf. Eq. Q. 



III. FITTING COEFFICIENTS 



As an illustrative example we study the T- fraction representation of Dawson's integral, which is closely related to 
the plasma dispersion function of a non-degenerate plasma, cf. Eq. §§§ and see also Ref. [12j ■ For a detailed discussion 
of Dawson's integral and its continued fraction representation see Refs. The application of T-fractions to this 

problem has already been considered by J. McCabe, see Ref. [VJ\ . We repeat a few of these results for illustration. 
Since Dawson's integral docs not have the proper expansion at z = oo, we consider the modified function 
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for a real variable x. This function does have the expansions 
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Using these expansions, one can calculate the T- fraction representation for foo{x) with the help of the modified Q-D 
algorithm, see App. [Cl In general, excellent convergence is found. The relative error as a function of the variable x 
for the nth approximant is shown in Fig. [TJ For n = 10, the error is better than 10 -6 for all x. 

As a second example, we discuss the highly degenerate plasma with 9 — > 0. Again, g(x) as given by Eq. ([8]) has to 
be supplemented by an imaginary part to ensure the correct asymptotic behavior 
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For this function, the expansions read for x — > oo and x - 
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(19) 



(20) 



respectively. Again, performing the modified Q-D algorithm leads to partial numerators and denominators as listed 
in Tab. IIIII As before, the overall convergence is good, which is verified by inspection of the relative error between 
the continued fraction representation and the analytic expression Eq. (fl7|) , see Fig.[2j Only in the vicinity of x = 1, 
the agreement drops to about 2% for the 10th approximant of the continued fraction. 

Next, we turn to the general case of arbitrary degeneracy. Typical expansion coefficients Hi are compiled in Tab.[l] 
Also, the expansion coefficients for large x are given essentially by the Fermi integrals for order 1/2. They are listed 
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e 


C-3 


C-5 


C-7 


C-9 


0.1 


1.387(-1) 


6.368(-2) 


3.897(-2) 


2.810(-2) 


0.2 


1.537(-1) 


8.275(-2) 


6.256(-2) 


5.856(-2) 


0.3 


1.749(-1) 


1.130(-1) 


1.071 (-1) 


1.305(0) 


0.5 


2.270(-l) 


2.837(-l) 


5.277(-l) 


1.228(0) 


0.8 


3.153(-1) 


4.418(-1) 


8.723(-l) 


2.509(0) 


1.0 


3.771(-1) 


6.042(-l) 


1.565(0) 


5.580(0) 


1.5 


5.359(-l) 


1.250(0) 


4.777(0) 


2.532(1) 


2.0 


6.979(-l) 


2.143(0) 


1.845(1) 


7.639(1) 



TABLE II: Expansion coefficients for the large x expansion of the function g(x), cf. Eq. ([5j). 

in Tab. [TT] for various values of the degeneracy parameter 9. Finally, the imaginary contribution due to the modified 
expression in Eq. (fT"3]) has to be determined. It reads 

i x 2 e"! x 4 e^ (e 7 * - 1) x 6 



lmg(x) w ^(l + e 17 ) - 
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(21) 



Using these coefficients, the T-fraction approximation is determined using the modified Q-D algorithm described 
in App. [Cj The resulting partial numerators and denominators are recorded in Tab. IIII1 The general form of the 
T-fraction considered here is 



Re 



Ho 



| x — ib\ 

with i being the imaginary unit, hq = 2/3. 
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FIG. 2: The relative error for a T-fraction approximation of the function fo(x), see Eq. (|17p . as a function of the variable x for 
the nth approximant. 



e 







0.1 


0.2 


0.3 


0.5 


0.8 


1.0 


1.5 


2.0 


bx 


4.244( 


-1) 


4.280(-l) 


4.393(-l) 


4.571(-1) 


5.023(-l) 


5.749(-l) 


6.217(-1) 


7.307(-l) 


8.283(-l) 




4.244( 


-1) 


4.280(-l) 


4.393(-l) 


4.571(-1) 


5.023(-l) 


5.749(-l) 


6.217(-1) 


7.307(-l) 


8.283(-l) 


&2 


9.234( 


-1) 


9.407(-l) 


9.900(-l) 


1.061(0) 


1.224(0) 


1.460(0) 


1.608(0) 


1.938(0) 


2.223(0) 


a 3 


4.522( 


-1) 


4.545(-l) 


4.652(-l) 


4.872(-l) 


5.458(-l) 


6.375(-l) 


6.981(-1) 


8.374(-l) 


9.592(-l) 


b 3 


9.613( 


-1) 


9.932(-l) 


1.071(0) 


1.184(0) 


1.430(0) 


1.731(0) 


1.953(0) 


2.413(0) 


2.782(0) 


0,4 


4.702( 


-1) 


4.730(-l) 


4.792(-l) 


5.074(-l) 


5.879(-l) 


6.695(-l) 


7.673(-l) 


9.678(-l) 


1.117(0) 


b 4 


9.790( 


-1) 


1.046(0) 


1.135(0) 


1.286(0) 


1.650(0) 


1.634(0) 


2.168(0) 


2.956(0) 


3.388(0) 


Or, 


4.804( 


-1) 


5.072(-l) 


4.800(-l) 


5.116(-1) 


6.525(-l) 


2.700(-l) 


6.921(-1) 


1.239(0) 


1.409(0) 


h 


9.878( 


-1) 


1.223(0) 


1.177(0) 


1.260(0) 


2.046(0) 


3.122(-1) 


1.545(0) 


3.720(0) 


5.740(0) 




4.865( 


-1) 


7.067(-l) 


4.471(-1) 


3.718(-1) 


9.483(-l) 


-4.755(0) 


7.390(-l) 


2.021(0) 


3.815(0) 


h 


9.923( 


-1) 


3.713(0) 


1.108(0) 


5.128(-1) 


5.835(0) 


-4.275(0) 


-8.132(-1) 


8.054(-l) 


-2.160(0) 


(17 


4.902( 


-1) 


3.404(0) 


-2.049(-l) 


-8.853(-l) 


5.056(1) 


8.397(-l) 


4.991(0) 


1.492(-2) 


-2.648(0) 


l>7 


9.949( 


-1) 


-4.507(-l) 


4.077(-l) 


-6.900(-l) 


-8.928(-l) 


9.416(-1) 


5.566(0) 


7.319(-3) 


2.244(0) 


ag 


4.926( 


-1) 


4.519(-1) 


2.686(0) 


7.086(-l) 


-8.812(-1) 


-2.490(0) 


-1.334(-3) 


7.783(1) 


1.698(0) 


bs 


9.965( 


-1) 


1.041(2) 


-2.871(0) 


7.651(-1) 


7.850(0) 


-2.190(0) 


-1.653(-3) 


7.817(1) 


-1.833(0) 



TABLE III: Partial numerators a* and denominators 6, for the T-fraction approximation of g(x) as a function of the degeneracy 
parameter 8, see Eq. (f22jl . Note, that bi — a-z. 



IV. CONCLUSIONS 

It has been shown that T-fractions are a powerful method to obtain reliable and compact approximative expressions 
for the dielectric response of ideal plasmas at arbitrary degeneracy. In particular, a set of eight partial numerators 
and denominators are sufficient to obtain a satisfying accuracy for most practical applications. Fast implementations 
for the dynamic collision frequency and the self-energy in GW approximations are possible based on these results 
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due to a notable acceleration compared to a direct evaluation of the integral representation. An extension to bosonic 
systems appears to be straightforward and is work in progress. 
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APPENDIX A: SMALL x EXPANSION COEFFICIENTS 

In order to regularize the hyper-singular integrals in Eq. ([4]), we consider an expansion of the Fermi function around 
zero 

1 1 c" y 2 _ (e" - 1) y A _ e" (l - 4e" + e 2 ") y e 

exp(y 2 /6-r]) + l ~ 1 + C -" (i + e'?) 2 T 2(l + e") 3 # 2 ~ " 6(l + e f ') 4 ^ 

e , '(-l + ll?/-lle 2 '' + e 31 ') y s 

24(l + e'0 5 ¥ + (A1) 

The calculation has been done using the computer algebra software MATHEMATICA [l8[ . In detail, we obtain for 
the integrands 

i i - c y 2 / 9 

h 2 (y) = — , (A2) 

zyy ' l + en 1 + e« / e ~v 

x (1 + e-") (l - e^ 9 ) + (l + C y^ 9 -^ y 2 /6 

h M = yt^ (l + c-") (i + c^-„) ' (A3) 

AKy! 2 (1 + e") 3 (e" + ev 2 / e ) L \ v l v l ) 

-e y2 ' 9 (2 - 2y 2 /9 + y 4 /6 2 ) + c" +y2 / e (-4 + 2y 2 /0 + y 4 /6 2 ) + c 2 " (2 + y 2 /9 + y 4 /6 2 )} , (A4) 

h 5 (y) = r? r [-6 (1 + c") 3 (-1 + e y2/e ) 

{ ' 6(l + e") (e" + ev 2 / e ) V y 1 \ ) 

+6 (1 + c r ') 2 (c" + c' y2/e ) y 2 /0 + 3 (-1 + c 2r ') (e" + c y2/9 ^j y 4 /6 2 + (l - 4e" + e 2?7 ) j/ 6 /*? 3 ] , (A5) 

with Hi = dy hi(y)/y 2l ~ 2 . Also, for evaluating the integrals involved, it is convenient to determine the expansion 
of each integrand for small x by reading off the corresponding terms in Eq. (|A1[) . Note, that the limits of Hi in the 
non-degenerate and highly degenerate case are analytically known from the corresponding expansions of Eq. and 
Eq. ijHJ), respectively. Note also, that the values of Hi and H 2 are given by Arista and Brandt while all higher 
expressions are given here for the first time. 

APPENDIX B: LARGE x EXPANSION COEFFICIENTS 

The Fermi integrals of order 1/2 involved in the large x expansion have been determined by MATHEMATICA [l8[ 
exploiting the relation of the Fermi integral to the polylogarithmic function [l3T ] 

F a (x) = -r(l + a)Li 1+Q (-e*) . (Bl) 

The Fermi integrals F_i/ 2 , F1/2 j ^3/2 > F5/2 can be easily obtained with the interpolation formulas given by Antia [l9l ]. 
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APPENDIX C: MODIFIED Q-D ALGORITHM 



We summarize the modified Q-D algorithm introduced by de Andrade et al. [15[ . We start with the series expansions 
Eq. (fTTj) and Eq. (TP2"]) and allow for zero coefficients /z,. Next, an auxiliary expansion is introduced by subtracting the 
expansion of the function Kj (1 + z) with a free parameter K from the above given expansion generating a new set of 
coefficients 7$. For a convenient choice of K, all 7, will be non-zero and the ordinary Q-D algorithm can be applied. 
The continued fraction representation 



fJ-o 



z-(5\ 



q< 0) z 



„(°) . 



W 



(CI) 



having the original series expansions Eq. (jllj) and Eq. (|12j) can be reconstructed from the coefficients 7^ by the 
following algorithm: 



• Determine 7^ from a given set of fii and a convenient choice of K as 

7i = M: 

• Perform the Q-D algorithm 



{-\) l+1 K £ 



with the initialization 

and i = 2, 3, ...,?' = ... , 
• Determine for i = 1,2,... 

&i+l,0 

with the initialization 



Ai r) = , £}' 



7r 
7r-l 



(C2) 

(C3) 
(C4) 



-2,-1,0,1,2,.... 

the following quantities leading to and 
- ^-2-(Af) + 4 0) ) , 



= UiX 



(0) 



-ai,o<5 i+ i , 

u i+ i + Wi - di-i,i_2 - 1) - Uj - 6»-i,»-2 ~ 1) 

-1 



1 



di,0 
Wi 



Wi+l 



1 



Q'l = , (3l = 



(0) 



7o - A-<5i 
ioi = ^0 , "i = 7o , 

60,-1 = , d -i = 



(o) ' 



(C5) 



APPENDIX D: NUMERICS OF CONTINUED FRACTIONS 



There are several algorithms to numerically calculate the value of a continued fraction. Here, we use the backward 
recurrence algorithm which is known for its numerical stability 0. Given n partial numerators a, and partial 
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denominators hi, we determine an approximation to to the continued fraction by setting t n to a small, non-zero 
number and perform the iteration 



u-i = jrrr , (di) 

bi + U 



for i = n, . . . , 1. 
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